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An  anode-supported  micro-tubular  solid  oxide  fuel  cell  (SOFC)  is  analyzed  by  a  two-dimensional  axisym- 
metric  numerical  model,  which  is  validated  with  the  experimental  I-V data.  The  temperature  distribution 
generated  by  the  thermo-electrochemical  model  is  used  to  calculate  the  thermal  stress  field  in  the  tubular 
SOFC.  The  results  indicate  that  the  current  transport  in  the  anode  is  the  same  at  every  investigated  posi¬ 
tion.  The  stress  of  the  micro-tubular  cell  occurs  mainly  because  of  the  residual  stress  due  to  the  mismatch 
between  the  coefficients  of  thermal  expansion  of  the  materials  of  the  membrane  electrode  assembly.  The 
micro-tubular  cell  can  operate  safely,  but  if  there  is  an  interfacial  defect  or  a  high  enough  tensile  stress 
applied  at  the  electrolyte,  a  failure  can  arise. 

©  2009  Elsevier  B.V.  All  rights  reserved. 


1.  Introduction 

The  solid  oxide  fuel  cell  (SOFC)  is  one  of  the  most  promising  can¬ 
didates  for  power  generation  because  of  their  high  efficiency,  clean 
generation  of  electric  power  and  silent  work  [1,2].  Currently,  the 
SOFC  prototypes  have  proven  its  ability  to  achieve  high  electrical 
efficiency,  especially  for  a  hybrid  system  with  gas  turbine.  Several 
different  cell  designs  are  being  investigated  [1-6].  Compared  with 
the  planar  SOFC,  the  tubular  one  has  several  advantages  as  men¬ 
tioned  in  literature  [2,4-6].  Recently,  because  of  a  faster  start-up 
and  higher  volume  power  density,  the  anode-supported  micro¬ 
tubular  SOFC  has  been  rapidly  developed  [2,4,5].  The  present  main 
concerns  to  SOFC  commercialization  are  the  high  manufacturing 
costs  and  short  lifetime.  The  structural  failure  is  one  of  major  issues 
for  the  high  degradation  rates. 

Many  investigations  have  been  carried  out  on  the  thermo¬ 
electrochemical  modeling  of  the  tubular  and  planar  SOFC  [2-4, 7-9], 
however,  only  a  few  have  been  carried  on  thermal  stress  up.  Two 
types  of  thermal  stress  are  present  in  a  SOFC.  The  first  one  is  the 
residual  stress  due  to  the  mismatch  between  the  coefficients  of 
thermal  expansion  (CTE)  of  the  materials  of  the  membrane  elec¬ 
trode  assembly  (MEA).  Finite  element  method  has  been  employed 
to  simulate  the  cell  stress  [3,6,9,10].  The  second  one  is  thermal  mis¬ 
match  stress  due  to  temperature  gradient  at  operating  conditions. 
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Selimovic  et  al.  [3]  and  Nakajo  et  al.  [6]  have  coupled  their  thermo¬ 
electrochemical  model  with  a  finite-element  tool  to  simulate  the 
stress  in  the  cell  during  operation.  In  order  to  quantify  and  predict 
the  magnitude  of  the  thermal  stress  in  a  micro-tubular  SOFC,  a  more 
specific  model  should  be  developed. 

In  this  study,  a  two-dimensional  (2D)  axisymmetric  numeri¬ 
cal  model  is  set  up  to  analyze  an  anode-supported  micro-tubular 
SOFC.  The  model  presented  in  this  paper  covers  momentum,  elec¬ 
trochemical,  heat-species  transport  and  the  thermal  stress  analysis. 
The  simulation  of  the  I-V  characteristics  versus  the  experimental 
data  from  the  cell  performance  tests  is  carried  out.  The  main  pur¬ 
pose  of  the  present  work  is  to  study  the  characteristics  of  the  stress 
field  including  the  residual  stress  and  thermal  mismatch  stress.  The 
approach  is  to  couple  a  thermo-electrochemical  model  for  the  com¬ 
putation  of  the  temperature  field.  After  that,  the  temperature  field 
is  imported  to  the  structural  mechanics  module  and  the  stress  field 
is  solved.  Temperature  gradients  and  thermal  stress  between  the 
components  are  evaluated.  The  Weibull  approach  is  developed  to 
investigate  the  effects  of  the  stress  on  the  MEA  reliability. 

2.  Experimental 

2.1.  Fabrication  of  anode-supported  micro-tubular  SOFC 

A  green  NiO/YSZ  tubular  support  is  prepared  by  an  extrusion 
method.  An  equal  weight  of  commercial  nickel  oxide  (NiO,  J.  T.  Baker 
Corporation)  and  8  mol%  yttria  stabilized  zirconia  (YSZ,  Tosoh  Cor¬ 
poration)  are  used.  After  the  extruded  tube  is  dried,  a  thin  film  of 
the  YSZ  electrolyte  is  coated  on  the  tube.  The  electrolyte-anode 
bilayered  tube  is  then  co-fired  at  1723  K  for  4  h.  A  Lao.8oSro.2oMn03 
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Nomenclature 

c  concentration  [ mol  m-3  ] 

Cp  specific  heat  capacity  [J  mol-1  K-1  ] 

Dy  diffusion  coefficient  [m2  s-1  ] 

E  Young’s  modulus  [GPa] 

F  Faraday’s  constant  [96485  C  mol-1  ] 

i  local  current  density  [A  nrr 2  ] 

i0  exchange  current  density  [A  nrr2  ] 

k  thermal  conductivity  [W  nr1  K-1  ] 

I<  permeability  [m2] 

m  Weibull  modulus 

ne  electron  transferred  per  reacting 

P  pressure  [Pa] 

Pf  failure  probability 

q  heat  source  [ W  nr 2  ] 

Q  volumetric  heating  rate  [ W  m-3  ] 

Rg  gas  constant  [8.314J  mol-1 I<-1  ] 

T  temperature  [K] 

Tre f  reference  temperature  [K] 

u  velocity  [ms-1] 

V0  reference  volume  [mm3] 

Greek  letters 

a  thermal  expansion  coefficient  [K-1  ] 

P  transfer  coefficient 

£rad  thermal  radiative  emissivity 

st  the  total  strain 

sth  thermal  strain 

r\  overpotential  [V] 

/a  dynamic  viscosity  [m2  s-1  ] 

p  density  [kg  m-3] 

a  conductivity  [S  m-1  ] 

cr0  Stefan-Boltzmann  constant  [5.67  x  10-8  W m2  K4] 
ex s o  characteristic  strength  [MPa] 

os  stress  [MPa] 

v  Poisson’s  ratio 

0  potential  [V] 

Subscripts  and  superscripts 
a  anode 

act  activation 

c  cathode 

ele  electrolyte 

j  cell  component  (anode,  electrolyte  and  cathode) 

rad  radiation 

rev  reversible 

t  total 

0  standard 


(LSM)/YSZ  composite  cathode  is  slurry-coated  on  the  YSZ  film,  and 
then  fired  at  1523 1<  for  2  h.  The  dimensions  of  the  cell  are  listed  in 
Table  1. 

2.2.  Single  cells  measurements 


Table  1 

Input  parameters  to  model  [2,11,12]. 


Descriptions 

Symbol 

Value 

Furnace  interior  diameter  (mm) 

6 

Tubular  interior  diameter  (mm) 

1.47 

Anode  thickness  (mm) 

0.8 

Cathode  thickness  (|jim) 

50 

Electrolyte  thickness  (p,m) 

10 

Anode  length  (mm) 

49.6 

Cathode  length  (mm) 

5.1 

Electrolyte  length  (mm) 

26.1 

Anode  conductivity  (S  irr 1 ) 

95  x  106T_1 

Cathode  conductivity  (S  irr 1 ) 

crc 

exp(-1150/T) 
42  x  106T_1 

Electrolyte  conductivity  (S  m- 1 ) 

CTele 

exp(- 1200/7) 
3.34  x  104 

Anode  thermal  conductivity  (W  m-1  K-1 ) 

ka 

exp(- 10300/7) 
1.86 

Electrolyte  thermal  conductivity  (W  m- 1 I<-1 ) 

ke  le 

2.16 

Cathode  thermal  conductivity  (W  irr 1  K-1 ) 

fee 

5.84 

Anode  density  (kg  m-3 ) 

Pa 

3310 

Electrolyte  density  (kg  m-3) 

Pe  le 

5160 

Cathode  density  (kg  irr3 ) 

Pc 

3030 

Anode  specific  heat  capacity  (J  kg-1 K_1 ) 

CP,a 

450 

Electrolyte  specific  heat  capacity  (J  kg-1 K1 ) 

Lp.ele 

430 

Cathode  specific  heat  capacity  (J  kg-1  K-1 ) 

Cp,c 

470 

Transfer  coefficient 

P 

0.5 

Electron  transferred  per  reacting 

ne 

1 

Permeability  (m2) 

K 

1  x  10-12 

Molar  entropy  change  of  cathode  reaction 

A  Sc 

81.63 

(J  mol”1  K"1 ) 

Thermal  radiative  emissivity 

£rad 

0.3 

the  cathode  side  is  connected  to  the  whole  cathodic  surface.  The  Ag 
wire  with  a  silver  mesh  is  used  for  current  collection.  The  perfor¬ 
mance  of  the  micro-tubular  cell  is  measured  after  the  open  circuit 
voltage  (OCV)  of  the  cell  reached  1.0  V. 

3.  Numerical  model 

The  model  makes  use  of  the  differential  equations  which  are 
integrated  in  the  2D  domains,  through  a  commercial  software, 
COMSOL  MULTIPHYSICS®.  Fig.  2  depicts  the  2D  schematic  of  the 
cell  used  for  the  numerical  simulation.  The  input  parameters  of 
the  model  are  shown  in  Table  1.  The  outputs  of  the  model  are  the 
distributions  of  temperature,  current  density  and  thermal  stress. 

3.1.  Model  assumptions 

•  The  model  is  based  on  steady  state,  using  humidified  H2  as  the 
fuel. 

•  Gas  flow  in  the  electrode  channel  is  a  fully  developed  laminar 
flow. 

•  The  reactant  gas  mixtures  are  approximated  as  ideal  gas  and 
incompressible.  The  fuel  cell  operates  with  100%  current  effi¬ 
ciency. 

•  The  electrochemical  reactions  are  considered  to  only  take  place 
at  the  electrode/electrolyte  interface. 

•  The  effect  of  gas  on  thermal  conductivities  of  porous  electrodes 
is  neglected. 


Fig.  1  shows  the  diagram  of  the  micro-tubular  cell  reactor.  The 
tubular  cell  is  located  at  the  center  of  a  furnace  with  a  furnace 
wall  temperature  of  1073  K.  The  inlet  fuel  is  composed  of  H2  and 
H20  (0.97:0.03  in  molar  fractions),  flowing  at  a  volume  flow  rate 
of  20  cm3  min-1  at  room  temperature.  The  fuel  gas  is  preheated  to 
1073  K  before  entering  the  cell.  The  cathode  side  is  open  to  air.  The 
current  collectors  at  the  anode  side  are  attached  to  the  outer  sur¬ 
faces  of  both  ends  of  the  anode  tube  and  the  current  collector  at 


Fig.  1.  The  diagram  of  the  micro-tubular  cell  reactor. 
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Fig.  2.  Schematic  of  a  micro-tubular  geometry  in  the  axisymmetric  2D  model. 


3.2.  Electrochemical  model 

The  Butler-Volmer  kinetic  equations  are  generally  adopted  for 
expressing  the  local  current  density  distribution. 


ia  =  *’a;b  •  (Ph2)°'25(Ph2 


x0.5 


exp  p 


2ne  •  F  •  (0rev  —  10a  —  0ele I )  A 

) 


electrode. 

p(u  •  Vu)  -  /x( V2u)  +  VP  =  0  (4) 

VP  =  ^V2u  +  (j£)  u  (5) 

The  steady  state  diffusion  and  convection  equation,  Eq.  (6)  is  solved 
to  obtain  the  species  distribution  in  all  domains. 

V(-Dir  Vq  +  qu)  =  0  (6) 

In  this  model,  the  conductive  and  convective  equation,  Eq.  (7)  is 
used  to  calculate  the  temperature  distribution  in  all  subdomains: 

V(-k  VT)  +  p  CP  u  VT  =  (l  (7) 

Eq.  (8)  is  used  to  express  the  radiative  heat  transfer  between  the 
outer  surface  of  the  cell  and  the  furnace  wall. 


qrad=ff0-erad(r4-7'4al|)  (8) 

The  ohmic  heating  at  the  electrode  and  electrolyte  is  calculated  by 
Eq.(9). 

Q  =  or  ■  V20  (9) 


The  heat  of  the  electrochemical  reactions  entropy  change  is  esti¬ 
mated  by  Eq.  (10). 


9 rev,  electrode  — 


T  •  ASeiectrode  •  ^electrode 
2  F 


(10) 


The  ASeiectrode  denotes  the  molar  entropy  change  of  the  electrode 
reaction  at  the  local  reacting  temperature  and  pressure  [11].  The 
heat  of  the  ion  transfer  across  the  electrode/electrolyte  interfaces 
is  calculated  by  Eq.  (11),  where  qac t  can  be  found  from  Eq.  (12). 


-  exp 


-(1  -  P) 


Tle  •  F  •  (0rev  —  10a  —  0elel) 


(1) 


Qirr,  electrode  —  Pact,  electrode  ‘  ^electrode 
Pact  =  0rev  —  I0electrode  —  0elel 


(11) 

(12) 


(Po2/Po2)1/4 

1C  -  Vo  • 


1  +(Po2/Po2)1/2 


exp  p 


]He  •  F  ■  (0rev  —  10c  —  06 1 e I )  \ 


Rg-T 


-  exp 


-(1-/6) 


ne  ■  F  ■  (0 rev  —  10c  —  0elel) 


(2a) 


Po2(Pa)  =  4.9  x  103  exp  ^  (2b) 

In  Eq.  (1)  it  is  adopted  that  the  reaction  order  is  0.25  for  H2 
and  0.5  for  H20  [8].  The  cathodic  exchange  current  density  can  be 
formulated  from  the  Ref.  [7],  where  p ^  can  be  found.  Since  some 
of  the  parameters  could  not  be  obtained  directly  from  experiments 
or  published  literature,  fitting  processes  have  to  be  involved  in  this 
study.  In  Eqs.  (1 )  and  (2a),  i^20  and  i° 2  are  the  adjustable  parameters 
to  fit  the  experimental  data.  i^20  is  7500 Am-2  for  the  anode  and 
i°20  is  2200  Am-2  for  the  cathode.  0rev  can  be  obtained  from  Nernst 
equation.  Electron  and  ion  charge  conservation  at  the  electrode  and 
electrolyte  is  expressed  by  Eq.  (3). 

— V(cr .  V0)  =  0  (3) 

3.3.  Thermo-fluid  equations 

The  steady-state  incompressible  Navier-Stokes  equation,  Eq.  (4) 
is  solved  to  calculate  the  velocity  distribution  in  the  gas  chan¬ 
nel,  and  the  Brinkman  equation,  Eq.  (5)  is  used  for  the  porous 


3.4.  Structural  mechanics  model 

Thermal  stresses  in  the  three-layered  structure  (anode,  elec¬ 
trolyte  and  cathode)  are  calculated  using  the  structural  mechanics 
module  of  the  commercial  finite  element  software  COMSOL.  The 
temperature  field  can  be  imported  into  the  structural  mechanics 
module  and  the  stress  field  can  be  solved.  In  Table  2,  the  Young’s 
modulus,  the  thermal  expansion  coefficient  and  the  Poisson’s  ratio 
of  the  cell  materials  used  in  the  study  are  given,  which  are  assumed 
to  be  constant  over  a  range  of 273-1723  K.  All  materials  are  assumed 
to  be  isotropic. 

The  value  of  the  zero  stress  temperature  is  very  important  for  the 
determination  of  the  magnitude  of  the  stress.  Many  uncertainties 
remain  in  this  area.  A  thin  film  of  the  YSZ  electrolyte  is  deposited 
on  the  anode  tube  at  a  temperature  of  1723  K.  At  this  tempera¬ 
ture  both  the  anode  and  the  electrolyte  film  are  stress-free  in  the 
present  work.  A  LSM/YSZ  composite  cathode  is  fired  on  the  YSZ  film 
at  1523  K.  At  this  temperature  the  cathode  is  stress-free. 

In  this  model  when  hydrogen  is  introduced  and  NiO  is  reduced 
into  Ni,  no  shrinkage  is  considered.  That  has  currently  been  reported 
in  open  literatures  during  NiO/YSZ  reduction  [10,18].  The  CTEs 
of  NiO/YSZ  and  Ni/YSZ  are  considered  to  be  equal  in  this  paper. 
This  was  also  considered  in  open  literature  [10].  When  the  above- 
mentioned  two  assumptions  are  adopted  in  this  model,  it  is 
acceptable  to  only  calculate  the  stress  distribution  for  the  Ni/YSZ 
anode  in  the  model. 

Finally,  the  thermal  stress  equation  is  solved  at  the  cell  oper¬ 
ating  temperature.  The  fixing  of  the  tubular  SOFC  is  modeled  by 
constraining  the  axial  displacement  at  the  fuel  inlet  boundary.  The 


Table  2 

Material  properties  of  the  cell  components. 
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other  boundaries  are  free  and  the  analysis  is  fully  elastic.  Ther¬ 
mal  loads  are  introduced  into  the  model  according  to  the  following 
equations,  Eqs.  (13)  and  (14). 


&th  —  &  (T  Tre f) 


(13) 


<?s  =  D  (£t  —  £th) 


(14) 


According  to  the  generalized  Hook’s  law,  the  following  elasticity 
matrix  D,  Eq.  ( 15)  is  used  for  an  isotropic  material  at  the  2D  axisym- 
metrical  model. 


D  = 


E 

(1  +  v)(l  -  2v) 


-  v  v  v 

v  1  -  v  v 

V  V  1  -  V 

0  0  0 


0 

0 

0 

1  -2v 
2 


(15) 


3.5.  Failure  probability  calculation 


The  rupture  often  occurs  at  the  MEA  due  to  the  brittle  nature 
of  the  ceramic.  For  this  reason,  it  is  valuable  to  analyze  the  fail¬ 
ure  probability  of  the  MEA.  The  Weibull  approach  is  considered  to 
estimate  the  risk  of  rupture  of  the  MEA  in  this  model.  The  Weibull 
analysis  [6]  calculates  the  failure  probability  (Pf)  of  the  MEA  with 
the  tensile  stress: 

Pfj  =  1  -  n-P  (-  /  C0*/  0s,o)mdVj/Vjo\  (16) 

where  crso,  m  and  V0  represent  the  characteristic  strength,  the 
Weibull  modulus  and  the  reference  volume,  respectively.  The  dif¬ 
ferent  materials  have  a  different  Weibull  modulus,  which  can  be 
obtained  experimentally.  A  greater  Weibull  modulus  indicates  that 
the  material  has  higher  reliability.  The  Weibull  modulus  of  the 
anode  used  in  this  study  was  adopted  of  Ref.  [6].  The  subscript  i  (1, 
2,  and  3)  denotes  three  principal  stresses,  o\,  cr2  and  or3.  The  sub¬ 
script  j  denotes  the  cell  components,  including  anode,  electrolyte 
and  cathode.  The  input  Weibull  parameters  for  SOFC  materials  used 
in  this  model  are  listed  in  Table  2. 


Fig.  3.  Comparison  of  simulated  and  experimental  performance  data. 

above-mentioned  three  r  dimensions.  Fig.  4  shows  the  r  direction 
component  of  current  density  along  fuel  flow  direction  at  the  differ¬ 
ent  r  dimension  (r  =  0.935  mm,  1.135  mm,  1.335  mm)  at  0.7  V.  From 
Fig.  4,  it  is  found  that  the  r  direction  components  of  current  den¬ 
sity  close  to  current  collectors  are  a  large  negative  value  at  all 
three  dimensions.  The  above  phenomena  show  that  the  current 
mainly  flows  from  the  outer  surface  nearby  the  current  collectors 
to  the  anode  inside,  and  the  flowing  direction  points  to  the  central 
symmetrical  axis.  The  absolute  value  of  current  density  gradually 
decreases  from  1.335  mm  to  0.935  mm.  These  results  indicate  that 
the  r  direction  components  of  the  current  density  close  to  the  cur¬ 
rent  collectors  are  bigger  than  those  at  other  positions.  After  the 
current  flows  into  the  anode,  the  r  direction  components  of  cur¬ 
rent  density  all  become  zero.  When  the  current  reaches  the  anodic 
reactive  zone  at  the  middle  of  cell,  the  current  value  becomes  pos¬ 
itive.  The  phenomena  show  that  the  r  direction  components  of  the 
current  density  flow  into  the  electrolyte. 


4.  Results  and  discussion 

4.1.  Model  validation 

The  cell  performances  obtained  from  the  model  are  compared 
with  the  experimental  performance  data  as  shown  in  Fig.  3.  The 
simulation  results  show  that  the  calculated  polarization  curves 
agree  very  well  with  the  experimental  data.  This  demonstrates  that 
the  present  model  is  accurate  enough  for  studying  the  behavior  of 
a  tubular  SOFC  operating  with  humidified  hydrogen. 

4.2.  Current  density  distribution  in  the  anode 

To  understand  the  current  flow  in  the  anode,  this  paper  gains 
three  current  densities  at  the  different  dimensions  which  are  par¬ 
allel  to  z  direction  and  are  evenly  distributed  in  the  anode.  Their 
distances  to  the  central  symmetrical  axis  are  0.935  mm,  1.135  mm, 
1.335  mm,  respectively.  To  compare  the  current  flow  at  the  differ¬ 
ent  anodic  locations,  the  model  calculates  the  r  and  z  direction 
component  of  current  density  along  fuel  flow  direction  at  the 
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Fig.  4.  The  r  direction  component  of  current  density  along  fuel  flow  direction  at  the 
different  r  dimension  (r = 0.935  mm,  1.135  mm,  1.335  mm)  at  0.7  V. 
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Fig.  5.  The  z  direction  component  of  current  density  along  fuel  flow  direction  at  the 
different  r  dimension  (r = 0.935  mm,  1.135  mm,  1.335  mm)  at  0.7  V. 

Fig.  5  shows  the  z  direction  component  of  current  density 
along  fuel  flow  direction  at  the  different  r  dimension  (r = 0.935  mm, 
1.135  mm,  1.335  mm)  at  0.7  V.  From  Fig.  5,  it  is  found  that  the  z 
direction  components  of  current  density  begin  to  appear  close  to 
the  current  collectors  at  all  three  different  dimensions.  Because  the 
current  flow  direction  at  fuel  outlet  is  opposite  to  that  at  fuel  inlet, 
the  current  density  value  at  outlet  is  negative.  After  the  current 
flows  into  the  anode,  the  z  direction  current  components  at  both 
ends  rapidly  increase  and  the  absolute  values  at  all  three  different 
dimensions  are  almost  equal.  From  Figs.  4  and  5,  in  the  process  of 
transferring  current  to  the  reactive  zone  of  the  cell,  the  r  direction 
components  of  current  density  are  zero  and  the  absolute  values  ofz 
direction  current  density  are  the  same.  It  indicates  that  the  current 
transport  is  even  at  every  investigated  position.  When  the  current 
reaches  the  reactive  zone  of  the  cell,  the  absolute  values  of  the  z 
direction  current  sharply  decrease.  Compared  with  Fig.  4,  it  can  be 
seen  that  the  decrease  is  because  the  z  direction  current  compo¬ 
nents  are  transferred  to  the  r  direction  current  components  and 
then  flow  into  the  electrolyte. 

4.3.  Temperature  distribution  in  the  cell 

Fig.  6  shows  the  temperature  change  curve  along  anode  outer 
surface  including  anode/electrolyte  interface  under  the  different 
load  conditions.  From  Fig.  6,  it  is  shown  that  when  the  cell  voltage 


Fig.  6.  The  temperature  change  curve  along  the  anodic  outer  surface  including  the 
anode/electrolyte  interface  at  the  different  load  condition. 
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decreases  from  0.9  V  to  0.5  V,  the  temperatures  at  the  reactive  zone 
of  the  cell  gradually  increase  and  the  temperature  maximum  is  at 
the  cell  center  along  z  direction.  The  maximum  reaches  1093  K  at 
0.5  V.  This  result  indicates  that  although  the  cell  locates  at  a  con¬ 
stant  temperature  furnace,  because  plentiful  heat  is  released  at  the 
electrochemical  reaction,  the  cell  temperature  is  always  higher  than 
the  furnace  temperature  at  the  cell  operating  condition.  The  cell 
temperature  rapidly  decreases  at  the  inert  zone  of  the  cell. 

Fig.  7  shows  temperature  gradients  distributions  in  the  cell  at 
0.7  V.  From  Fig.  7,  it  is  shown  that  firstly,  the  temperature  gradient 
gradually  increases  from  both  ends  of  anode  to  the  cell  center,  and 
the  maximum  temperature  gradient  is  at  the  double  ends  (points  B 
and  C)  of  the  reactive  zone  of  the  cell,  and  that  secondly  the  tem¬ 
perature  gradient  gradually  decreases  from  the  double  ends  of  the 
reactive  zone  of  the  cell  to  the  middle  of  the  cell.  From  Figs.  6  and  7, 
it  is  found  that  the  temperature  gradient  at  the  position  of  the  tem¬ 
perature  maximum  is  smaller  than  at  other  positions  of  the  reactive 
zone  of  the  cell. 

4.4.  Thermal  stress  analysis 

Fig.  8  shows  the  arrow  plot  of  the  principal  stress  at  0.7  V.  The 
arrows  express  the  principal  stress  directions,  and  the  arrow  sizes 
are  proportional  to  the  principal  stress  value.  The  part  enclosed  by 
the  circle  in  Fig.  8(a)  is  scaled  up  to  Fig.  8(b)  and  the  part  enclosed 
by  the  circle  in  Fig.  8(b)  is  scaled  up  to  Fig.  8(c).  The  inward  arrow 
expresses  the  compressive  stress;  the  outward  represents  the  ten¬ 
sile  stress.  In  Fig.  8,  it  is  shown  that  the  principal  stress  is  much 
larger  in  the  electrolyte  than  in  other  parts.  The  principal  stress  in 
the  electrolyte  is  mainly  compressive  along  the  z-axis. 

The  cell  presents  some  material  and  geometrical  singularities 
at  its  edges  A,  B,  C,  and  D  in  Fig.  2  according  to  the  Refs.  [10]  and 
[14].  Because  of  the  high  stress  level  close  to  the  above-mentioned 
edges,  these  positions  in  the  cell  are  potentially  harmful  for  the 
cell  structure.  Fig.  9  shows  the  first  principal  stress  distribution  at 
0.7  V.  In  strength  theory,  principal  stress  is  stress  in  the  principal 
plane  which  is  free  of  shear  stress.  There  are  three  principal  stresses, 
cr \ ,  a 2  and  a3.  By  comparing  the  value,  the  order  is  determined, 
cti  >  (T2  >  cr3-  In  Fig.  9,  the  stress  in  the  anode  is  tensile  below  7  MPa 
except  at  the  cell  singular  points.  The  stresses  in  the  cathode  and 
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Fig.  8.  The  arrow  plot  of  the  principal  stress  at  0.7  V.  The  part  enclosed  by  the  circle 
(a)  is  scaled  up  to  (b),  and  then  the  part  enclosed  by  the  circle  in  (b)  is  scaled  up  to 
(c). 

electrolyte  are  compressive,  primarily  due  to  the  lower  CTE  of  the 
cathode  and  the  electrolyte  materials.  The  maximum  of  the  first 
principal  stress  locates  close  to  the  both  ends  of  the  electrolyte, 
points  A  and  D,  due  to  material  and  geometrical  singularities.  From 
the  Ref.  [10],  the  cell  singularities  are  harmless  under  the  proper 
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cell  edge  structure,  so  the  singularities’  effect  on  the  cell  failure  is 
not  considered  in  this  paper. 

Figs.  10  and  11  show  the  second  and  third  principal  stress  at  0.7  V. 
In  Figs.  10  and  11 ,  it  is  shown  that  the  stress  is  strongly  compressive 
about  270  MPa  in  the  whole  electrolyte.  The  compressive  stress  in 
the  cathode  is  less  than  10  MPa.  It  is  known  that  compared  with 
compressive  stress,  tensile  stress  is  easier  to  cause  ceramic  material 
fracture.  The  above  results  indicate  that  the  cathode  is  relatively 
stable  under  a  low  level  of  compressive  stress.  The  compressive 
stress  in  the  electrolyte  is  much  higher,  but  still  below  the  failure 
level  because  the  fracture  strength  of  electrolyte  at  the  compressive 
condition  has  been  reported  to  be  about  1  GPa  [13]. 

4.5.  Failure  probability  analysis 

The  specific  values  of  the  stress  field  cannot  provide  a  failure 
probability  at  any  particular  operating  point.  In  order  to  investigate 
the  effects  of  the  stress  on  the  MEA  reliability,  the  Weibull  approach 
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Fig.  11.  The  third  principal  stress  distribution  (Pa)  at  0.7  V. 
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Fig.  12.  Failure  probability  of  the  anode  and  MEA  with  various  voltage  of  the  cell. 

is  developed,  and  by  which  the  influences  of  the  cell  voltage  on  fail¬ 
ure  probability  are  researched.  The  influences  of  the  cell  voltage  on 
the  failure  probabilities  of  the  anode,  MEA,  cathode  and  electrolyte 
are  shown  in  Figs.  12  and  13,  respectively.  The  results  indicate  that 
failure  probabilities  of  the  all  parts  except  the  cathode  increase  with 
the  increase  of  the  cell  voltage.  The  probability  of  failure  is  lower  at 
low  load  in  the  cell  components  because  the  effect  of  the  temper¬ 
ature  is  stronger  than  that  of  the  temperature  gradient  in  the  cell. 
In  the  cathode,  on  the  other  hand,  the  probability  of  failure  is  lower 
at  high  load  when  the  effect  of  thermal  gradients  is  stronger  than 
that  of  the  temperature.  However,  it  has  to  be  mentioned  that  the 
results  in  this  model  are  limited  to  a  particular  choice  of  material 
properties. 

Weibull  analysis  only  calculates  a  positive  value.  With  regard 
to  the  calculation  of  the  stress  field,  the  tensile  stress  is  positive 
and  the  compressive  is  negative.  The  stress  in  the  anode  is  mostly 
tensile;  the  stresses  in  the  cathode  and  the  electrolyte  are  mostly 
compressive.  So,  the  results  of  the  probability  of  a  failure  of  the 
cathode  and  the  electrolyte  are  much  smaller  than  that  of  the  anode. 
From  Fig.  12,  it  is  found  that  the  failure  probabilities  of  the  anode  and 
MEA  are  at  the  same  order  of  magnitude.  So,  it  can  be  considered 
that  the  failure  probability  of  the  anode  is  the  decisive  factor  of  the 
failure  probability  of  the  MEA.  As  shown  in  Figs.  12  and  13,  the  dif¬ 
ference  of  the  failure  probabilities  between  the  cell  operating  state 
and  the  open  circuit  is  less  than  3%.  From  Figs.  7  and  9,  it  can  be  also 
found  that  when  the  maximum  temperature  gradient  has  attained 
8000 Km-1,  the  anodic  tensile  stress  is  still  below  7 MPa.  These 
phenomena  indicate  that  the  micro-tubular  cell  can  resist  the  high 


Fig.  13.  Failure  probability  of  the  electrolyte  and  cathode  with  various  voltage  of  the 
cell. 


Fig.  14.  The  cell  strain  energy  density  (J  rrr3 )  at  0.7  V. 

thermal  shock  and  that  the  cause  of  stress  is  mainly  the  residual 
stress  due  to  the  mismatch  between  the  CTE  of  the  materials  of  MEA. 

The  residual  stress  is  caused  by  the  different  CTEs  of  the  elec¬ 
trolyte  membrane  and  electrodes.  When  the  CTEs  of  all  components 
are  the  same  and  the  temperature  is  constant,  no  stress  can  be 
found  in  the  cell.  In  the  anode-supported  SOFC,  the  CTE  difference 
between  the  anode  and  electrolyte  membrane  is  the  decisive  factor 
of  the  stress  in  the  cell.  Compared  with  the  electrolyte,  the  CTE  of  the 
anode  is  easy  to  adjust.  Because  the  CTE  of  the  anode  is  higher  than 
that  of  the  electrolyte,  an  effective  method  to  decease  the  failure  of 
probability  of  the  cell  is  to  decrease  the  CTE  of  anode. 

Though  the  above  analysis  has  proved  that  the  tubular  cell  can 
operate  safely  enough,  the  Weibull  approach  is  only  applied  at  ten¬ 
sile  stress  condition.  Failure  at  interfaces  between  electrolyte  film 
and  anode  substrate  will  happen  if  the  film  is  subjected  to  a  suffi¬ 
cient  high  compressive  stress  so  that  the  elastic  strain  energy  stored 
in  the  film  exceeds  the  energy  required.  Actually,  the  high  com¬ 
pressive  stress  can  cause  the  layer  to  fail  when  an  interfacial  defect 
occurs  or  a  critical  tensile  stress  is  applied  on  the  layer.  In  this  work, 
the  risk  of  failure  has  been  investigated  by  calculating  the  elastic 
strain  energy  stored  in  the  compressive  thin  layer.  Fig.  14  shows 
the  cell  strain  energy  density  at  0.7  V.  The  strain  energy  stored  in 
the  cathode  layer  is  only  about  44  kj  m-3 ,  whereas  the  energy  stored 
in  the  electrolyte  layer  reaches  1575  kj  m-3.  This  value  exceeds  the 
energy  requirement  of  the  electrolyte  film  (ca.  1000  kj  m-3  [  10]).  So, 
the  failure  can  arise  if  there  is  an  interfacial  flaw  or  a  high  enough 
tensile  stress  applied  at  the  electrolyte. 

Several  limitations  of  this  simulation  method  exist.  Such  as,  in 
the  SOFC,  the  electrodes  are  porous.  If  a  micro-scale  electrode  can 
be  added  in  this  macroscopic  model  and  a  multiscale  model  is  used, 
the  model  can  show  more  information  and  more  accuracy  about  the 
stress  distribution.  At  the  process  of  the  anodic  reduction,  several 
assumptions  are  considered  in  this  model.  After  the  reduction,  the 
micro-morphology  of  anode  is  very  important  for  the  cell  perfor¬ 
mance.  The  micro-morphology  change  from  NiO/YSZ  to  Ni/YSZ  can 
be  taken  into  account  in  detail  in  the  future. 

5.  Conclusions 

A  2D  SOFC  model  is  constructed  in  this  study  to  evaluate  the 
thermal-fluid  behavior,  electrochemical  reaction  and  thermal  stress 
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in  a  SOFC  single  cell.  The  temperature  distribution  is  the  inputs  for 
a  finite  element  software,  which  is  used  to  estimate  the  thermal 
stress  field.  The  experimental  data  are  employed  to  validate  the 
simulation  results.  Good  agreements  have  been  observed  between 
simulation  results  and  experimental  data.  The  current  flow  phe¬ 
nomena  in  the  anode  are  analyzed  in  detail.  The  result  indicates 
that  the  current  transport  in  the  anode  is  the  same  at  every  inves¬ 
tigated  position.  The  temperature  gradient  at  the  position  of  the 
temperature  maximum  is  smaller  than  at  other  positions  of  the 
reactive  zone  of  the  cell.  The  micro-tubular  SOFC  can  resist  the  high 
thermal  shock  and  the  stress  occurs  mainly  because  of  the  residual 
stress  due  to  the  mismatch  between  the  CTE  of  the  materials  of 
MEA.  It  is  not  recommended  to  keep  a  quite  slow  rate  of  sintering 
temperature  at  the  process  of  the  SOFC  manufacture,  especially  at 
high  temperature  because  the  stress  is  smaller  close  to  the  stress- 
free  temperature.  The  effective  method  to  decease  the  failure  of 
probability  of  the  cell  is  to  decrease  the  CTE  of  anode  for  anode 
supported  SOFC.  Matching  the  CTEs  is  the  most  important  factor 
to  attain  a  higher  reliability  of  the  cell.  Through  the  Weibull  anal¬ 
ysis  it  has  been  proved  that  the  tubular  cell  can  operate  safely,  the 
cell  failure  can  arise  if  there  is  an  interfacial  flaw  or  a  high  enough 
tensile  stress  applied  at  the  electrolyte.  In  the  process  of  SOFC  oper¬ 
ation,  one  should  pay  attention  to  impose  the  external  force  as  low 
as  possible. 
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